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ABSTRACT 

The solar convective zone, or SCZ, is nearly adiabatic and marginally con- 
vectively unstable. But the SCZ is also in a state of differential rotation, and its 
dynamical stability properties are those of a weakly magnetized gas. This ren- 
ders it far more prone to rapidly growing rotational baroclinic instabilities than 
a hydrodynamical System would be. Thèse instabihties should be treated on the 
same footing as convective instabihtes. If isentropic and isorotational surfaces 
coïncide in the SCZ, the gas is marginally (un)stable to both convective and rota- 
tional disturbances. This is a plausible resolution for the instabilities associated 
with thèse more gênerai rotating convective Systems. This motivâtes an analysis 
of the thermal wind équation in which isentropes and isorotational surfaces are 
identical. The characteristics of this partial differential équation correspond to 
isorotation contours, and their form may be deduced even without précise knowl- 
edge of how the entropy and rotation are functionally related. Although the exact 
solution of the global SCZ problem in principle requires this knowledge, even the 
simplest models produce striking results in broad agreement with helioseismology 
data. This includes horizontal (i.e. quasi-spherical) isorotational contours at the 
pôles, axial contours at the equator, and approximately radial contours at mid- 
latitudes. The theory does not apply directly to the tachocline, where a simple 
thermal wind balance is not expected to be valid. The work presented here is 
subject to tests of self-consistency, among them the prédiction that there should 
be good agreement between isentropes and isorotational contours in sufïiciently 
well-resolved large scale numerical MHD simulations. 
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1. Introduction 

The problcm of diffcrcntial rotation in the solar convective zone (SCZ) has vexed the- 
orists for many years. At the simplest level, one is faced with the fact that the only région 
of the sun characterized by significant turbulence is also the only région showing any signif- 
icant differential rotation. This is a cogent reminder of the hazards of modeling a rotating, 
turbulent gas by an effective viscosity parameter. The SCZ is far too complex to yield to 
such an approach. Indeed, large scale numerical simulations have been in place for many 
years, trying to elicit the rotation profile of the sun's outer layers (Brummell, Cattaneo, & 
Toomre 1995; Thompson et al. 2003). While thèse studies are now at a stage where they 
can begin to make contact with the observational data, they have not yet been able to re- 
produce the salient features of the observed solar rotation profile. This is especially true for 
the isorotational contours (hereafter "isotachs"), which tend to be cylindrical in simulations, 
but are decidedly noncylindrical in the sun. 

In this paper, we examine an cffect that could be important for understanding the 
SCZ rotation profile, but has not been sufficicntly cmphasized in previous investigations: 
the extrême sensitivity of ionized differentially rotating gas to the présence of even very 
weak magnetic fields of arbitrary geometry (Balbus & Hawley 1994, Ogilvie 2007). Far from 
behaving as a passive vector field, a weak magnetic field triggers rapidly growing unstable 
local modes in rotating Systems that would be hydrodynamically unstable. The field endows 
the gas with degrees of freedom that have no hydrodynamical counterpart. The classical 
manifestation of this behavior takes the form of what is known as the magnetorotational 
instability (MRl), widely regarded as the process responsiblc for turbulence and enhanced 
transport in accrction disks (Balbus & Hawley 1991, 1998). But destabilization by a weak 
magnetic field is a process that extends beyond the domain of rotationally supported disks. 
Such fields are agents of dynamical destabihzation in baroclinic and convective flows as well 
(Balbus 1995), with potentially interesting applications to the SCZ. 

It is csscntial to understand two key countcrintuitivc points at the outsct. The first is 
that destabilization by a magnetic field can occur even when the field plays essentially no 
rôle in the dynamical equilibrium: a weakly magnetized gas does not behave "almost hydro- 
dynamically" with respect to its stability properties. The second is that the destabilization 
is independent of the field strength and field geometry. In MRI simulations, the disruption 
remains vigorous and self-sustaining even in highly convoluted turbulent flow. This is what 
makes MHD processes potentially so important for understanding rotating, stratified flows. 

It is the emphasis on thèse points that sets the présent approach apart from earlier 
work that also addressed the stabihty of rotating, stratifled, magnetized flows (e.g. Acheson 
1983, Ogden & Fearn 1995). In thèse earlier studies, the focus is upon toroidal flelds. 
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nonaxisymmetric disturbances, and instability driven directly by the magnetic field. By 
contrast, what is crucial here are poloidal fields and axisymmetric disturbances. Moreover, 
while the instability of intcrcst relies on the présence of a magnetic field, the (nonconvective) 
seat of free energy is the differential rotation, not the field itself. The rôle of the field is to 
provide the critical degrees of freedom to the fluid required to tap into the (destabilizing) 
differential rotation. 

We shall begin our présentation, however, not with magnetic fields, but with a theo- 
retical observation that may be taken at a purely phenomenological level, independently of 
deeper MHD considérations. This is the remarkable agreement with the overall pattern of 
solar isorotation contours produced by a simple calculation in which a dominant thermal 
wind balance from the vorticity équation is combined with the assumption that entropy and 
angular velocity gradients are counteraligned. Whatever the cause of the counteralignement 
may be, the analysis on its own is highly suggestive. (A connection between entropy and 
angular velocity gradients has in fact been advocated on purely hydrodynamical grounds 
[Miesch, Brun, & Toomre 2006].) This material is presented in §2. In §3, we discuss in 
some détail the case for an MHD coupling between angular velocity and entropy gradients. 
The final section is a critical discussion of unresolved issues provoked by this work, and a 
summary. 

2. Solving the thermal wind équation 
2.1. Review 

Let {R, (p, z) be a standard cylindrical coordinate System, and (r, 6, (j)) a standard spher- 
ical coordinate System. Unit vectors are denoted Br, eg, etc. The angular velocity Q, is 
assumed to be independent of 0, but otherwise gênerai. Our notation for the fluid variables 
is likewise standard: v is the velocity, P is the gas pressure, p is the mass density, and B is 
the magnetic field. 

We adopt a fiducial value of 2.5 x 10~^ rad s~^ for the angular velocity Q, in the SCZ, 
corresponding to rotation velocities between 1 and 2 km s"^. This is well in excess of the 
~ 30 m s~^ expected of convective velocities, but such a direct comparison is not necessarily 
the most relevant one. A more telling comparison is between the squared Brunt-Vàisàlà 
frequency 

(we have adopted a value of 10~® for d\YïPp~^ /dlnr [Schwarzschild 1958]) and the rotational 



1 dPdlnPp 

p'-f dr dr 



3.8 X 10"^^ s-2 



(1) 
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parameter 



"^^^ 1.8 X 10-^2 (2) 



d\nR 

This is evidently about 5 times as large as A^^. (Whether R or z is used in the angular 
velocity gradient is not critical at midlatitudes.) Taking this estimate at face value, the 
rotational parameter is significantly larger than N'^. The significance of this will shortly 
become apparent. 

Let us first considcr Systems whosc cquilibrium velocity is differential rotation, v = 
Rfl{R, z)e(i,. The thermal wind équation foUows directly from the time-steady form of the 
vorticity équation, 

^^ = ;^(VPxVp).e^. (3) 
Expressing the right side in r, 9, spherical coordinates: 

_ 1 fdpdP _dpdP\ 

Now, rewrite this in terms of the entropy gradients: 

ao^ _ 1 rdPds dPds\ 

~d7 ~ Cppr \ dë'd^ 'd^'dë)' 



where S is the spécifie entropy, 

k 



- In Pp + constant, (6) 



k is the Boltzmann constant, and Cp is the constant pressure spécifie heat. 



c.= [:^^)k (7) 

In the SCZ, the r gradient of P clearly dominâtes over the 9 gradient, whereas the r gradient 
of S is unhkely to greatly exceed the 9 gradient. (In fact, we shall presently argue just the 
opposite.) We may then conclude that 

^ _^dS ^ g djlnPp-^) 
dz Cpr de -fr de ^ ^ 

where g = —{1/ p)dP/dr is the gravitational field in hydrostatic equilibrium, ignoring the 
small effects of rotation. Significant latitudinal entropy gradients are required to avoid 
cylindrical isotachs. 
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Equation ([8]) is known as the thermal wind équation (Kitchatinov & Riidiger 1995, 
Thompson et al. 2003). This équation holds in our MHD analysis, because we explicitly 
assume that the magnetic fields are sufïiciently weak in the equilibrium state that they do 
not affect the large scale rotation profile. (For a 10 G field and a density of 0.05 g cm^^, 
the Alfvén velocity is 13 cm per second.) In full spherical coordinates, the thermal wind 
équation is 



Even before a detailed solution is developed, a simple scaling argument reveals something 
of interest here. The r gradient of S is generally determined by the need to transport the solar 
luminosity by thermal convection (Schwarzschild 1958, Clayton 1983). On the other hand, 
the 6 gradient is, from équation ([8]), linked directly to differential rotation. If the fiducial 
numbers ([T]) and ([2]) are reasonably accurate, at midlatitudes équation ([8]) implies that the 6 
gradient of S will significantly exceed the r gradient. The point of interest here is that this 
anisotropic feature is directly seen in the Q gradient deduced from the helioseismology data. 
This suggests that there is a deeper dynamical coupling présent between S and fl, beyond 
just the gênerai trend that one goes up as the other goes down. We will develop this idea 
more fuUy in the next two sections. 



The thermal wind équation is a familiar tool to practitioners of solar rotation theory. 
The gross features of the sun's angular velocity profile {Q decreasing polewards from the 
equator) can be understood relatively simply with the aid of this équation and some reason- 
able assumptions of the efïiciency of convection in the présence of Coriolis forces (Thompson 
et al. 2003). The idea is that convection in the equatorial direction is impeded by Coriolis 
forces, resulting in a more efficient transport of heat along the rotation axis. The pôles are 
then régions of higher spécifie entropy compared with the equatorial zone, and the resulting 
6 gradient in 5* drives axial gradients in Q. 

This is a well-motivated and physically plausible scénario. It seems natural therefore, 
to search for explicit model solutions to the thermal wind équation incorporating this idea. 
Moving upward from the equator, we expect the angular velocity to decrease as the entropy 
increases. Clearly the gradients of thèse quantities are in broadly opposite sensés. But 
the final paragraph of §2.11 suggests that "counteraligned" may be a better description then 
"broadly opposite." This raises the question of what the solutions to équation (jH]) would 
look like if the two gradients were in fact precisely oppositely aligned. 




(9) 



2.2. Analysis 
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FoUowing the notion that isentropic and isorotational surfaces coincide, we assume that 
S = S'(fi^). (5* should not change if fl changes sign, hence the dependence on Q'^.) Equation 
then becomes 

+ = 0' (10) 



dr \C pr'^ sin 6 cos 6 r J d6 
where S' = dS/dVt^. The solution of équation (fTUI) is that Vt^ is constant along the charac- 
teristic 

d6 ta.n6 qS' , , 

— = (11) 

dr r Cpr^ sin^ cos0 

But if Ç)? is constant along this characteristic, then so must be 5*'. This is therefore a self- 

contained, ordinary differential équation for ^ as a function of r, precisely the isorotational 

contours that we seek. To solve équation (fTTj) . let y = sm9. Then, our differential équation 

simplifies to 

_ _2gS^ (12) 
dr r Cpr'^ 



Multiplying by and regrouping, 

d{r^y^) 2gS' 2GMqS' 



(13) 



dr Cp Cpr"^ 

where G is Newton's constant and Mq is a solar mass. (We have ignored the local self-gravity 
of the SCZ.) Since 5" is a constant along the characteristic, this intégrâtes immediately to 

sin^ e = R^ = A-- (14) 
r 



where A is a constant of intégration and 



2GMqS' 



B = (15) 

We have inserted a minus sign since 5* will generally be a decreasing function of Q^. Indeed, 
as will become very clear, the characteristics make no sensé if S' is positive, but a very great 
deal of sensé if it is négative. In this model, the solar isotachs are given by a remarkably 
simple formula. 



To estimate the magnitude of B, note that it may be written 
B (2GMQ/rQ\ ( d\ïiPp-^'\ ( dhir 



(16) 



y 7r|i72 J \ d\nr J ydlnVl"^ 

The first factor is large, of order 10^. The second, we have already estimated at 10^^. The 
helioseismology data suggest that the third factor ranges between 1 and 10, and is larger 
near the equator. Crudely speaking, we expect B/r^ to be of order unity or less. There 
should not be a large différence in scale between A/r^ and B/rQ. 
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2.3. Alternative isotach solution 

By way of comparison, consider the solution of the thermal wind équation that would 
obtain were the angular momentum counteraligned with the entropy, rather than the angular 
velocity. (In hydrodynamic baroclinic turbulence, one might expect the coupling to be of 
this nature since entropy and angular momentum tend to be retained by a displaced fiuid 
élément.) Then S = S{P), where / = R^Q is the spécifie angular momentum. We dénote 
dS/dP' by Si2. Our analysis proceeds along lines identical to §2.21 and P satisfies the PDE 

The contours of constant are found to be of the form 

1 1 . 5« 

m = = A + — 18 

where Ai is an intégration constant, and now 

B, ^ (19) 

Note that TqBi is a dimensionless constant of order unity. We will use this solution as a 
point of comparison in the next section. 



2.4. An explicit solution 

Let us turn to the data to see how our solutions fare, beginning with our first, Q based 
solution. As an illustrative example, consider the reduced problem in which 5" is constant, 
not just along a particular characteristic, but everywhere. Then B is also everywhere con- 
stant. If Q is now specified on some particular radial shell as a function of 6, we may write 
down the solution everywhere. For this purpose, it is easiest to choose the solar surface 
radius r = Tq. 

Let the angular velocity at r = r© be Qq{cos^ 6q). We use the fit of Ulrich et al. (1988): 

fi0(cos2eo) = 27r (451.5 -65.3 cos^ ^0-66.7 cos^^o) nHz. (20) 

Note that 9q carries a subscript to indicate that it is the particular value of 9 along each 
trajectory characteristic (ÎT^ that intersects the surface shell r = Tq. Equation (fl^ becomes 

sin^ ^ = sin^ ^0 + 5 - , (21) 
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or 

where 

w = r/rQ^ x^ = w'^sm^6. (23) 

(Equations fl2T]) and fl22|) are also valid if i? is a function of 0086*0.) Substituting équation 
fl22|) for cos^ ^0 iiito équation fl2U|) then générâtes the solution everywhere in the shell. 

The resuit of this procédure is shown in figure (1) for the cases B/r^ = 0.5,0.6. The 
detailed fit of the contours to the actual data is not perfect — there is no tachocline in thèse 
simple models, and the true high latitude isotachs show stronger curvature, foUowing spher- 
ical shells, before turning upward. But the overall trend of the isotachs being prédominant ly 
quasi-spherical at high latitudes, increasingly radial at midlatitudes, and axial at small lat- 
itudes is unmistakable. Moreover, fitting our solution to the observed surface data, while 
convenient, is unlikely to show off its best form: thermal wind balance probably breaks down 
near the solar surface (Thompson et al. 2003). Given the simplicity of our direct approach, 
the qualitative agreement is both striking and encouraging. 

The iso-angular-momentum contours of équation (fTSjl can also be used to construct an 
explicit solution. In this case, the surface angular momentum / fit is 

/(cos^ 6*0) = 27rrQ sin^ ^0 (451.5 - 65.3 cos^ ^0 - 66.7 cos^ ^0) . (24) 

Instead of équation (l22l) . we have 

cos 6*0 = 1 — sin^ ^0 = 1 — 

Substitution of ( l25l) into (HM générâtes the full solution for the spécifie angular momentum 
l{r,6), and the angular velocity solution follows immediately from 

^(cos^ 9o) = {271 /x^) sin^ ^0 (451.5 - 65.3 cos^ 9o - 66.7 cos^ ^o) (26) 

In figure 2, we show two représentative diagrams of the isorotational contours taken 
from this alternative angular momentum based approach. In gênerai the contours are too 
cylindrical, to some extent exhibiting the same syndrome often seen in numerical SCZ sim- 
ulations. The contrast between figures 1 and 2 is very apparent. There seems to be a real 
linkage between 5* and Q, and it matters very much that the coupling is between S and 
Q, not 5* and l. It is possible that the refractory nature of the cylindrical contours of the 
simulations is due to an S" — / coupling that remains too strong, as noted in in §2.31 While 



+ reBi 1 



1 



—1 



(25) 
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Fig. 1. — Contour plot of Q{r,d) in solar interior, using surface fit of Ulricli et al. (1988). 
SCZ boundaries marked in white. Calculation based on eqs. [20] and characteristic équations 
of isorotational contours, [22]. B/r^ = 0.5 (left), 0.6 (right). 




Fig. 2. — Contour plot of Q{r, 9) in solar interior, based on counter aligned entropy and 
angular momentum gradients, using surface fit of Ulrich et al. (1988). SCZ boundaries 
marked in white. Calculation based on eq. [26] and characteristic équation of iso-angular- 
momentum contours, [Mj. TqBp = 0.2 (left), 0.5 (right). 
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it is possible that this may be cured by a more highly resolved treatment of the turbulent 
fluid, it is also possible that magnetic fields may be playing a non-negligible rôle, enforcing 
an S — Q coupling by field line tethering of the fluid éléments. We pursue this possibility in 
§3. 

In the remainder of the paper, we focus exclusively on our original S ~ Q solution. 



2.5. Tightening the contours 

By allowing the B parameter to vary from one isotach to another, the isorotational 
contours we have found can become more tightly spaced near the pôles. In this sensé, our 
solutions admit, but do not demand, something reminscent of tachoclinic structure. Equation 
(ÎT^ may be written 

r = ^/-^ (27) 

i-{RyA) ^^^^ 

If the variation of Q leads to a nearly constant ratio for B/A close to the outer radius of the 
radiative zone Vrad, but very différent A values, tachoclinic structure results. In the polar 
régions near the axis of rotation, the contours would ail converge to Vrad, while corresponding 
to very différent Q values. This would look very much like a tachocline. 

Does this make physical sensé? Equation fl?I]) implies A = r^sin^^o + -S/r©, or 

B _ r© 
A 1 + (r| sin^ Ôo/B) ' 



(28) 



assuming that the data are initially specified on r = r©. Hence, B oc sin^ 9o would produce 
tachoclinic structure. 

Physically , this would mean that B oc dS/ dVL^ is small and négative near the pôle, and 
that the entropy decreases toward the equator as Vl increases. In the vicinity of the equator 
the entropy drops sharply. This is not unreasonable behavior: near the pôle, unencumbered 
by Coriolis déviations, convection is most effective, whereas at the equator, the opposite is 
true. At the same time, the observations suggest that Vt is changing rapidly near the pôle, 
and only slowly at the equator. This is consistent with our simple picture. 

We return to équation (j22D and replace B /r^ by 

5/r| =r7i+r72sin2^o. (29) 



Mathematically, this corresponds to the first two terms in a Taylor séries expansion of S' 
as a function of siv? 6q (or equivalently cos^^o)- By varying r/i and rj2 we can go between 
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a singular tachocline rji = and the previous case Î72 = 0. Substituting fl29|) in fl22|) and 
solving for cos^ 6*0 gives 

cos 00 = 1 ■ — . (30) 

w + r]2{to - 1) 

Using équation (1301) in (1201) now produces the interior structure Q{r,9). 

With two parameters available, one can of course produce somewhat better fits to the 
rotation profile. In practice, the improvement over §2.41 is noticeable, but not dramatic. In 
figure 3, we show on the left the isotachs for rji = 0.3, 772 = 0.2. This gives a very respectable 
fit to the data away from the tachocline, r > O.TSr©, say. If we wish to include an explicit 
tachocline in our modeling, our earlier considérations suggest that we should restrict ourselves 
to small values of rji = 0. On the right side of figure 3, the interesting case of 77 = 0.12 and 
ri2 = 0.8 is presented. A striking "tachocline" structure appears, though formally it lies just 
beneath the SCZ. The solar midlatitude radial contours and equatorial cylinders régions are, 
however, rather well-represented. Note as well the gentle nonmonotonic behavior of Q{R) 
near the equator, a feature seen in the helioseismology data. Increasing the value of eta2 to 
bring the tachocline to larger radii (while keeping the surface layers fixed) appears to cause 
too much global distortion, though an exhaustive parameter search has not been performed. 

At this stage, thèse results are suggestive, but not more than that. It seems likely 
that more complex choices for S"(f2^) could improve the contour fits, but the agreement is 
impressive even in the simplest models. A truly compelling explanation of the tachocline 
will involve more than just the thermal wind équation and the outer convective zone layers. 
The solar tachocline arises from the complex coupling of the rigidly rotating radiative core 
and the overlying strongly shearing convective zone, not the demands of the surface rotation 
and vorticity conservation, and very différent dynamical processes are likely to be involved. 
Our simple, prescription valid above the tachocline (r ^ 0.75rQ) may represent a sort of 
outer SCZ solution that asymptotically matches onto an inner solution in which tachocline 
dynamics become locally dominant. 



2.6. Generic features of isotachs 

The isotachs we have found have a very distinctive "viking helmet" structure. This 
unusual feature is also characteristic of solar isorotation contours, and worth examining in 
isolation. 

To understand the gênerai structure of the isotachs, rewrite eqation (1271) as 

r B/Atq a 



Tq 1 - (r|/A)(/?2/r|) 1 - (3R^/rl 



(31) 
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which defines our two parameters a and 

a = B/rç,A, (3 = rlj A (32) 
We now transform to Cartesian coordinates in a méridional plane, 

y = {r /tq) cosO, x = (r/r©) sin 6^. (33) 
The Cartesian contour structure is slightly easier to visualize: 

1/2 

(34) 

We must of course view équation (l34l) through the convective zone slot 0.7 < \J x"^ + y"^ < 
1. There are three types of contours. The first is a typical polar contour, emerging perpen- 
dicular to the axis before bending upward, seen at high latitudes in figure 1. The second 
class, visible at midlatitudes in figure 1, is hidden in the radiative zone at small x, and does 
not émerge into the SCZ until the contour is well separated from the axis, at which point 
the isotach has a quasi-radial character. The third contour class is deeply buried, running 
along a very small radius in the core (not seen), completely disappearing in the bulk of the 
sun {y is imaginary). It then makes a sudden leap upward into the equatorial région of the 
SCZ, where it appears nearly axial. Thèse are the low latitude régions of figure 1, corre- 
sponding physically to Taylor columns of constant rotation on cylinders. Equation (l34ll thus 
displays the three key traits one associâtes with isotachs: horizontal near the pôles, radial 
at midlatitudes, and cylindrical near the equator. 

3. Axisymmetric modes in a rotating, baroclinic, weakly magnetized gas 

3.1. Preliminaries 

The starting point of the analysis of §2 was that the dominant balance of the vorticity 
équation is given by the thermal wind relation. We adopted a phenomenological connection 
between 5* and VL that succeded in reproducing the observed gênerai behavior of the solar 
isorotational contours. In this section we examine a possible reason for this coupling. We 
suggest that the underlying cause of the S—VL coupling is to be found in the gênerai dynamical 
stability properties of a magnetobaroclinic fiuid. 

We have already noted that in accretion disk applications, the combination of magnetic 
fields and differential rotation is extremely destabilizing, even if the field is very weak and 



a 



j3x 



2^2 
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not particularly well-ordered. This is a conséquence of the MRI. The analysis présentée! here 
is an extension of the accretion disk problem, but a highly significant extension: we con- 
sider the most gênerai dynamical axisymmetric response to a weakly magnetized barochnic 
System. To avoid possible confusion, we will reserve the tcrm "MRI" to apply only to the 
magnetic destabilization process in rotation-dominated disks. The topic of this section is 
the "magnetobaroclinic instability." 

The System we analyze is a proxy that shares important features with the sun. It 
consists of a body of self-gravitating gas that has arbitrary axisymmetric angular velocity 
and cntropy profiles. Instability in the form of turbulent thermal convcction is trcatcd on 
the same footing as rotational instability. Both the thermal and rotational profiles can 
in principle be altered by turbulent fluxes arising from magnetobaroclinic instability, since 
neither profile is fixed by the requirements of hydrostatic equilibrium. This differs from the 
behavior of an accretion disk, whose rotational profile, generally Keplerian, is not at liberty 
to change. 

It is an elementary fact that a convectively unstable stratified gas tends to alter its 
thermal gradient to a nearly adiabatic configuration, thereby regulating the linear instability 
itself. What is novel here is that we extend this notion to include simultaneously both the 
rotational and thermal rcsponses. This is generally not something investigators havc pursued, 
because at first sight it docs not appcar to be particularly promising. Hydrodynamically, the 
differcntial rotation of the sun is not close to instability. It is only when magnetic ficlds are 
considcrcd that rotational instabilitics arc raiscd to the same Icvcl as convcctivc instabilitics. 
Indeed, even a uniformly rotating magnetized gas is just "marginally stable," tipping to 
instability with only a slightly adverse angular velocity gradient. 

The invocation of magnetic instability may strike some readcrs as dubious. Is a turbulent 
convective zone fertile ground for process that rcquircs at least some dcgrcc of field cohérence? 
Is it justified in the analysis to prescribe Q a priori when in fact it is built up by convection? 
Thèse questions can ultimately be settled only by well-designed numerical simulations. In 
the meantime, let us first understand the behavior of our proxy magnetic System, a challenge 
in itself. We will then be in a better position to address more thorny issues. 

3.2. Stability of a magnetized baroclinic gas 

We seek to understand the linear stability properties of a gas in which convective, rota- 
tional, and magnetic effects are treated as co-equals. As before, let {R, 0, z) be a standard 
cylindrical coordinate System. Consider an axisymmetric rotating gaseous body whose equi- 
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librium angular velocity fi and entropy S are allowed to be functions of both R and z. 
We consider local perturbations of plane wave form, exp{ik • r — iujt). Here, k is the local 
axisymmetric wavenumber, r the position vector, uj the angular frequency (or growth rate) 
and t is the time. Such disturbances satisfy the dispersion relation (Balbus 1995, Balbus & 
Hawley 1998): 

fc2 



^4 I 2 



7P W 



-AVL\k-VA) = 



(35) 



where va is the Alfven velocity 



B 



Va 



and 



dz dR 



LU 



{k-VA?. 



(36) 



We follow the stability arguments of Balbus (1995). The variable cc^, and hence cj^, 
must be real. We may therefore détermine stability by noting those conditions under which 
uj"^ passes through zéro. The solution c<j^ = is possible if 

[k-VA? = ^ (^4fi2 + -l(PP)P(lnPp-^) + 1^V{R'Q')^ . (37) 



To assure stability, this équation cannot have any solutions for (fc-t?^)^, hence the right side 
must satisfy 

(38) 



4fi' + —(VP)V(lnPp-^) + ■^V(R''n^) < 0, 
7P il'* 



a condition that does not involve the magnetic field, though it pertains only to a magnetized 
fluid! Notice in particular that the field geometry is unimportant. The hydrodynamical 
stability condition, by way of contrast, would be 



— (VP)V(\nPp-^') + -^V(R^Q^) < (hydrodynamic stability), 
7P R'^ 

an altogether différent and far more easily satisfied requirement. 



(39) 



If, in équation fl38|) . we set x = k^/k^ and expand the T> operator, we may recast the 
inequality as 



x'^Nf + X 



dP d In Pp"^ ^dPdhi Pp-^ 



where 



7P \ dz dR dR 
1 dPdlnPp-^ 



dz 



R- 



,dœ_ 

dz 



dn^ 

d\nR 



> 0, 



n: 



P7 dz dz 



1 dPdlnPp- 

'^m dW 



(40) 



(41) 
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In the thermal wind models of §2, > 0, A^^ < 0, and N"^ + dfl'^/dlnR > throughout 
the bulk of the SCZ at midlatitudes. 

Two conditions will ensure that the quadratic poljTiomial in x is positive. We refer to 
thèse as the magnetized H0iland criteria, after the investigator who solved the corresponding 
hydrodynamic problem (Tassoul 1978). The first is 

since this means that either very large or very small x is positive. This criterion generally 
seems to be satisfied throughout the bulk of the convection zone. The second criterion foUows 
from requiring that the discriminant of the quadratic x polynomial fHUl) should be négative 
so that there are no real roots. The resuit of this somewhat lengthy calculation i^ 

~dl)\~dR 5i dR ^ 

In the course of deriving this resuit, we explicitly use the component of the vorticity 
équation, 

,dn'^ 1 fdpdP dpdP' 



^ dz p^\dRdz dzdRJ' ^^^^ 

which will be recognized as the starting point for our thermal wind analysis in §2. It can be 
shown that équations (H2!) and (143|1 together imply 

a somewhat more stringent requirement than équation ( l42l) by itself, but an obvions set of 
constraints that also follows simply upon inspection of équation fHUl) . 

The second magnetized H0iland criterion (H3!) states that the component of VS* X VTi 
should be positive in the northern hémisphère and négative in the southern hémisphère 
for stability. (See figure 4.) Marginal stability by this second important criterion, would 
correspond to entropy and angular velocity surfaces coinciding. In our problem, the gradients 
would be oppositely directed. It is to be noted that this criterion holds regardless of field 
geometry or strength, as long as the Alfvén velocity is relatively small. 

The marginalization of the linear instability must be viewed at présent as a plausible 
outcome of the induced turbulent flow, rather than a certainty. But assuming that the 



^This resuit was first derived using a variational approach by Papaloizou & Szuszkiewicz (1992). 
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Fig. 3. — As in figure 1, but calculation now based on équations [20] and [30] with rji = 0.3, 
r]2 = 0.2 (left), and rji = 0.12, 772 = 0.8 (right). The two-parameter fit on the left is a slight 
improvement over the earlier single parameter models. The fit on the right is striking in its 
overall resemblance to the SCZ, though the "tachochne" formally hes below the convective 
zone lower boundary. Equatorial and midlatitude contours are well-represented. 




Fig. 4.— An UNSTABLE (left) and a STABLE (right) alignment of and VQ for 
northern hémisphère disturbances. Marginal stability, used in our solution of the thermal 
wind équation ([9]) corresponds to précise counter-alignment of thèse gradients. 
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system ultimately does arrange its entropy and angular velocity gradients to curtail unstable 
magnetobarotropic modes, it would appear to do so by passing through gradually diminishing 
values of [k • va)"^- Eventually the wavelength will exceed the size of the physical domain, 
and there can be no question of an unstable mode at this point. In practice however, the 
relevant lengthscale is likely to be determined by the cohérence length of the magnetic field, 
not the size of the SCZ. This is probably a dimension not very différent from the largest 
convective eddy scale. 



3.3. Marginal magnetobaroclinic instability 

Section lÏÏ^ suggests that the underlying dynamical explanation of the near coïncidence of 
constant entropy and angular velocity surfaces, which seems to be a good phenomenological 
model of the helioseismology data, is marginal stability to axisymmetric magnetobaroclinic 
modes. This raises a number of questions. Why thèse modes? Why focus on axisymmetry? 
Might nonaxisymmetric modes be more unstable? 

The answer to the first question is that the modes considered are just the standard 
convective motions at the heart of solar turbulence, but which find themselves subject to 
magnetic fields and rotation. Marginal instability arguments and near adiabatic température 
profiles for stellar convection zones are uncontroversial. They are based on simple physical 
reasoning, not complex turbulence calculations. Here we suggest that there is an important 
augmentation to thèse arguments needed when both magnetic fields and rotation are présent. 
Because of the coupling introduced by the magnetic field, rotational instabihty must be 
considered with convection from the very start. 

To see why this might be so, as well as to answer the second and third questions, 
begin with the nonaxisymmetric dispersion relation for a magnetized, uniformly rotating 
gas. This may be derived using exactly the same procédure followed in Balbus (1995). For 
a wavenumber with (j) component m/R, the dispersion relation fISïïl) becomes 



2 



1 TD 



-4n''{k-VAf = (46) 



where 



uj' = u;'-{k.VAr, N' = Nl + N'^=-^^^^^^ (47) 

P'j or or 

and r is the spherical radius. Consider first a purely hydro dynamical rotating system, = 
0. For stability, 

N' > 0, {VP)V{\nPp-^) -An' <0 (48) 
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The first inequality, a requirement for nonaxisymmetric modes, is the standard Schwarzschild 
criterion. Ultimately a convectively System is characterized by a small négative value of A^^, 
just enough to maintain marginal levels of turbulence. The second inequality, required of ail 
axisymmetric modes, may be written 



This is guaranteed if 



1 fdPdlnPp-^ ^dPd In Pp-^ 



7P \ dz dR dR dz 



+ A^^ + 4f]2 > 0. (49) 



A^^ + 4fi2 > 0, Ni > (50) 



Thèse are slightly more restrictive than the nonaxisymmetric requirement, but hardly con- 
strain the rotation at ail. 

With the inclusion of even a weak magnetic field, there is a significant change in the 
stability of axisymmetric modes. The stability requirements of the dispersion relation (jl^D 
for any finite k • Va are 

A^^ > 0, Ni > 0. (51) 

The nonaxisymmetric requirement A^^ > is superfiuous. The présence of a magnetic field 
precludes any rotational stabilization. 

Finally, when we allow the combination of magnetic fields and differential rotation to 
be présent the axisymmetric modes are elevated to the rôle of key players. Now axisym- 
metric modes that would be stable by the Rayleigh criterion can be destabilized. Shearing 
nonaxisymmetric disturbances still couple to convective motion, but marginalizing the Brunt- 
Vâisâlâ growth rate is not enough to ensure dynamical stability. The rotation profile built 
up by convection must be marginalized as well (cf. équation [13]), and this is a much more 
stringent requirement in weak field MHD than in hydrodynamics. This is the central point 
of this section. 



4. Discussion and summary 

Our principle conclusion — that the rotational profile of the sun (and presumably other 
late type stars) is a magnetic phenomenon — is far reaching, and many readers may still be 
skeptical. Perhaps the most controversial points are the implicit assumption that the Alfvénic 
k • Va coupling remains vigorous on small scales in a turbulent fiuid, and that the gradients 
of entropy and angular velocity should be accorded equal respect in gauging the overall 
dynamical stability of the SCZ plasma. Thèse issues must ultimately be established or refuted 
by well-designed numerical MHD simulations. It bears emphasis, however, that there is a 



-19- 



precedent of nearly two décades of intensive numerical simulations of the magnetorotational 
instability in astrophysical disks. Global MHD simulations are vigorously unstable in the face 
of fully developcd turbulence and higlily convoluted fiow. The bulk of the magnetic energy 
is generally to bc found at large scales in a turbulent fluid (e.g. Fromang & Papaloizou 
2007), and the small wavelengths that drive instability are well-couplcd to the cmbedded 
fiuid magnetic field. The question is whether this is equally true of an MHD turbulent fiuid 
whose primary instability is convective, with the rotation profile built up by the convection 
itself. This important question remains open, at least for the moment. 

In any case, the rcsults of §3 of this papcr suggcst that the cntropy and angular vclocity 
in the SCZ know about each other one way or another. Bccause of the properties of ax- 
isymmetric instabilities in a weakly ionized gas, there is a dynamical basis for the belief that 
they may well be functionally related: important axisymmetric instabilities are controled 
when the entropy and angular velocity are counteraligned. But the fit of the thermal wind 
équation solution with the helioseismology angular velocity contours speaks for itself. It is, 
at the very least, empirical évidence for S ~ S'(n^) throughout much of the SCZ. 

Of course the fit is not completely perfect, nor should it be. The actual helioseismology 
data show polar contours hugging spherical shells before turning sharply outwards (Thomp- 
son et al. 2003). Our thermal wind contours are more smooth and Icss spherically curved. 
The data show contours closing ncar the surface at cquatorial latitudes, something that our 
characteristic-based theory docs not rcproduce. It is interesting as well to note that the the 
marginal stability arguments break down in this région, since dP/dz approachcs zéro and 
there is no longer a requirement of parallel entropy and angular velocity gradients. More- 
over, Thompson et al. (2003) have pointed out that thermal wind balance appears not to 
hold near in the outermost layers of the SCZ, where turbulent transport can no longer be 
neglected. In fact, our model results appear to be too smooth here, in comparison with the 
data. 

But our approach is very simple, and that it succeeds as well as it does is striking. 
If the magnetobaroclinic marginal stability arguments of this paper are not correct, then 
it is a coïncidence that such arguments lead to solutions of the thermal wind équation in 
broad agreement with the helioseismology data. If the arguments are correct, however, then 
hydrodynamical simulations are unlikely to reproduce the solar angular velocity contours, 
unless there is an as yct unknown purely hydrodynamical S—Q coupling. In the simulations of 
Miesch et al. (2006), when such a connection was put in "by hand" the contour fit noticeably 
improved. To reproduce the observed SCZ isotachs in numerical simulations without such 
forcing, the présence of a magnetic field may be essential, and care should be taken to ensure 
that the most rapidily growing magnetobaroclionic local instabilities are resolved. Internai 
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dynamics should take care of the rest: the System will evolve to counteract the instabilities 
to the extent it can. An excellent way for the flow to do this is to counteralign its entropy 
and angular velocity gradients. The resulting isotachs should then be in broad agreement 
with the helioseismology data. 

Capturing the local unstable modes may not be an easy task. To do so, wavenumbers 
k satisfying 

{kvAf ^ R-g^ (52) 

will need to be resolved (Balbus & Hawley 1998). This leads to X/tq ~ 10^'^ for typical SCZ 
values. Since a fraction of this wavelength must be represented on the grid (or higher order 
spectral basis functions), this is beyond the resolution of most of the simulations performed 
to date. Numericists may wish to consider artificially enhancing magnetic effects to drive 
the System toward MHD marginal stability as a possible means to improve the fit of the 
computed rotation contours. 

In principle, one may also test our claims by using the helioseismology data for Q{r, 9) 
and working backwards. From équation (fTUl) . the quantity 

sin e cos 9 [dn^/dr) [dn^/dO) - sin^ 9, (53) 

which is directly proportional to S"(fi^), should have the same isocontours as VL itself. But 
to test this would require a reliable extraction of the partial derivatives of Vt^. 

The search for an explanation of the SCZ isorotational contours has been long and not 
without frustration. If the work presented here is correct in its basic essentials, it would be 
a step forward. But if for some reason the consistency between the solutions of the thermal 
wind équation and the marginal stability requirements of local magnetobaroclinic modes is 
simply an accident, if weak field instabilities are ultimately not effective in the SCZ, even this 
réfutation would represent a form of progress by ruling out a viable alternative. Régulation 
of the SCZ turbulence by marginal stability to magnetobaroclinic modes is a well-posed, 
directly testable concept, involving a dynamical domain that is as yet underexplored. It 
merits serions considération. 
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